Reconstruction of Liouvillian Superoperators 
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We show how to determine (reconstruct) a master equation governing the time evolution of an 
open quantum system. We present a general algorithm for the reconstruction of the corresponding 
Liouvillian superoperators. Dynamics of a two-level atom in various environments is discussed in 
detail. 



Proper description of quantum dynamics of open sys- 
tems is essential for our understanding of physical pro- 
cesses in many areas of physics, starting from quan- 
tum optics to quantum cosmology. In general an open 
system can be represented as a system S interacting 
with an environment £ ||. In this paper we consider 
the archetypal system + environment model which is 
specified as follows: Let H s denotes a Hilbert space 
of the system S, and 7i E is the Hilbert space associ- 
ated with the environment E. The Hamiltonian H SE — 
H s ® 1 E + Hi n t + l s ® H E of the composite system S®E 
acts on TL S ® H E . It is assumed that S © E is a closed 
finite- dimensional system which evolves unitarily. The 
density operator p SE (t) of this composite system is gov- 
erned by the von Neumann equation with the formal solu- 
tion p SE (t) = cxp[-i(t-t )H SE }p SE (t Q )exp[i(t-t )H SE }, 
where the initial state is p 3E (to) = P s (to) © P E (to) and 
fi = 1. The reduced dynamics of the system S is then 
defined as 



Ps(t) : = T(t,t )p s {t Q ) = Tr E [p SE (t)} 



(1) 



By definition, T(t, to) is a linear map which transforms 
the input state p s (to) onto the output state p s (t). In this 
paper we address the question how to determine (recon- 
struct) the master equation which governs the time evo- 
lution of the reduced density operator p s (t). This master 
equation can be written in the convolutionless form Q 
(we omit the subscript S) 



dt 



p(t)=C(t,t )p(t)- 



(2) 



which is possible due to the fact that in the finite- 
dimensional Hilbert spaces matrix elements of density 
operators are analytic functions. Consequently, T{t,to) 
are non-singular operators (except may be for a set of 
isolated values of t) in which case the inverse operators 
T(t, to) exist and the Liouvillian superoperator can be 
expressed as 



£(t,t ) 



dt 



T(t,t ) 



t-\t,t ) 



(3) 



and 



We note that T(t,to) is uniquelly specified by H s 
by the initial state p E (to) of the environment. 

In this paper we propose a general algorithm how to re- 
construct the Liouvillian superoperator £(t,to) from the 



the knowledge of the time evolution of the system den- 
sity operator pit). In fact there are two aspects of this 
problem. Firstly, p(t) can be given as a result of a se- 
quence of quantum-tomography measurements || , such 
that at each time t the system density operator is re- 
constructed from the measured tomographic data. From 
these experimental data then the Liouvillian which gov- 
erns the open system can be reconstructed (see Exam- 
ple A). Secondly, the density operator pit) is determined 
from the knowledge of the unitary evolution of the com- 
posite S © E system [see Eq.(|l|)]. From this knowledge 
the master equation (^J) is determined. In both cases 
the dynamics of the open system is given exclusively in 
terms of the system operators. Environmental degrees of 
freedom are completely eliminated from the reduced dy- 
namics. Nevertheless, the state of the environment may 
change during the time evolution due to the interaction 
with the system. That is, we do not employ the assump- 
tion that the environment is a "big" reservoir which does 
not change under the action of the system (see Examples 
B and Cj. 

In order to reconstruct the Liouvillian superopera- 
tor C(t,to) we have to determine firstly the linear map 
T(t,to) given by Eq.(|l|). This part of the reconstruc- 
tion can be performed with the help of the algorithm 
recently proposed by Poyatos, Cirac and Zoller M. This 
algorithm works as follows: Let us assume that the 
system S has been initially prepared in a pure state 

|^(io)) = Sr=o c iil*i) wnere are basis vectors in 
the (N + l)-dimensional Hilbert space TL S of the sys- 
tem under consideration. It is further assumed that the 
environment is initially prepared in a state p E (to) — 
J2 ai a 2 d aia2 \a.i) E (a 2 \, where \a>i) E are basis vectors in 
the Hilbert space 7i E of the environment. 

In general, the physical process T(tk 1 to) is determined 
by a transformation acting on basis vectors of the sys- 
tem and the environment (in what follows we omit in all 
expressions the explicit reference to the initial time to) 



'i s ai 



T(t k 



N 



^2 X! E (n3l)(a 1 p 1 )(tk)\jl) s \/3l) 
31=0 /3i 



(4) 



The output density operator p(tk) of the system at time 
tfc is obtained when the transformation (|j) is applied to 
the initial state of the system-environment p(to)®p E (to) 
and then the tracing over the environment is performed, 
so that p(tk) can be written as 



1 



JV 



P(tk)= C »l( C *2)*%i/ i2 )(*fe): 



(5) 



%\ ,22—0 



where (JV + l) 2 operators R(i lt i 2 )(tk) are defined as 



A' 



R(i lt i 2 )(t k ) = £, (ii,<2)(ji,i2)(*fc)|jl)0'2|, (6) 

hda=0 

with 

■^(*i,»a)0'i Ja)(^k) = ^ Q i Q 2 

a l , a 2 ,7 

X %lil)(ai7)(*0^(i 2 j 2 )(a 27 )(*k)' ( 7 ) 

From Eq.(^|) it follows that the process l~{t k ) for a 
given time tk is completely determined by (JV + l) 2 op- 
erators R(i lt i 2 ){tk), which in turn are specified by the 
(JV + l) 2 x (JV + l) 2 matrix elements D( ilti2 )(j 1 lj2 )(t k ). 
We note that-R^^)^) have the properties 



(8) 



Tri?( ilii2 )(t fc ) = 5 iu i 2 ; 

("R(tl,t2 )(**)) = %2,il)(*fe)' 

or, equivalently, 

AT 

D (ii,i2)(ii,j 2 )(*fe) = DihMXhdi^k)- ( 9 ) 

We also note that neither R^ ili2 ^{t k ) nor -D(i 1 ,i 2 )(j 1 , : ; 2 )(ifc) 
depend on the initial state /3(to) 01 the system and for- 
mally they fulfill the conditions 



lim R (lul2) (t k ) = \h)(i 2 \; 
lim D {iui2)(jlj2) (t k ) = S iujl S i2j2 . 



(10) 



Poyatos et al. Q have shown that in order to spec- 
ify the (JV + l) 2 operators R/ ili3 \(tk) one has to con- 
sider (JV + l) 2 specific (see below) initial conditions 
= £^ =0 cj 1,fca) |i 1 ) where h, k 2 = 0, 1, JV 
and to measure the corresponding (JV + l) 2 output den- 
sity operators p^ kl ' k2 '(t k ) which can be expressed as 



JV 



R {klM) {t k ) = J2 ^ 1 ,fe2)(M,i2)Ail ) i2)(*fe), (11) 

ii 1*2=0 

where 

Mfrtox hAt) = # iM H<£ M) r. (is) 

If the (JV + l) 2 initial conditions |vE'^ 1,fe2 ^ ) ) ln are chosen 
so, that the matrix M^^iiiM) gi ven by Eq. ( |l2] ) is 



invertible, then the set of Eqs.(]ll|) can be solved with re- 
spect of the operators ■Rr» 1 ,i 2 )(ffc). Alternatively, one can 
express the matrix elements Du lt i 2 \rj 1 ,j 2 )(tfc) as functions 
of the in and out states of the measured system, i.e. 



X M (ii,i2)(k u k 2 ) S (k 1 ,k 2 ){ji,h)( t k), (13) 



iV 



fei,fes=0 



where the (JV + l) 2 x (JV+ l) 2 matrix S is defined as 

%, fc2 )(i 1 ,i2)(*fc) = 0'i|p (fcl ' fc2) (*fe)li 2 ). (14) 



The matrix M is the inverse of M and has the prop- 
erty Z)fc 1 ,fc 2 =0-^(ji,j2)(fci,fc2)-^(fci,fc2)(ii,i2) = ^'lji^jr 

So this is how the process T(t k ) can be reconstructed 
from the measured in and out states. To make the re- 
construction possible the matrix M has to be invertible. 
Obviously, there are many choices of such matrix. In 
particular, Poyatos et al. Q have proposed M given by 
Eq.(|i"2|) with complex amplitudes c^ 1 '^ specified as 



„(ki,ka) 



(8i,ki +<Ji,fe 2 )/v / 2 if hi > fc 2 
(£i,k! +^i,fc 2 )/v / 2 if fci < k 2 



(15) 



The reconstruction process described above gives us a 
set of operators Ru t> i 2 \(tk) which describe the transition 
of the system from the state p(to) to the state p(t k ) 
at a given time t k . In principle, one can perform a 
whole sequence of such reconstructions at different times 
ti,t%, ...Ik so that the reduced dynamics of the studied 
system can be reconstructed from the measured data. 

Now our task is to determine (reconstruct) from a set 
of measurements of the output states p( kl ' k2 '(t) for given 
input states p^ kl,k2 \to), the form of the Liouvillian su- 
peroperator C(t) in Eq.(||). To do so, we firstly note, 
that when the time evolution of the operators p( kl,k2 '(t) 
is governed by Eq.(^), then taking into account the ex- 
pression ([ll]) and the assumption that the matrix M is 
invertible, we find that the operators R^ lt i 2 ^(t) are also 
governed by the same master equation, i.e. 

j t R (lul2) {t) = C{t)R {lul2) {t), (16) 

with the initial conditions given by Eq.(nOf). Alterna- 
tively, taking into account the expression (ph we obtain 
from Eq. (Rq) a set of linear differential equations for ma- 



trix elements D 



(il,i 2 )(fel,fe 2 ) 



(*) 



d . . 

-^ D {iui2)(k u k 2 )\t) 



N 

E 

jl ,32=0 



D 



(»l.»2)(j'l,j2) 



32)(ki,k 2 ) 



(*), 



(17) 



with the initial conditions (|l0|). Here the matrix 



G 



Ul<j2)(ki,k 2 ) 



(t) is defined as 



2 



G 



(h,j2)(k 1 M)( t ) = ( fc il (^Wb'i)0'2 



(18) 



and it uniquely determines the Liouvillian superoperator 

We already know how to reconstruct matrices D from 
the measured data for arbitrary time t (from these data 
we can also evaluate the corresponding time derivatives). 
Providing the matrix Dr ilti2 \(j 1 j 2 \(t) is not singular its 

inverse Dq 1 j 2 \u 1 ^ 2 \(£) can be found and then the re- 
constructed matrix G^ 1 j 2 ^^ 1 ^{t) is given by a simple 



expression 



N 



G (jl,j2)(fcl,fc2)(*) 



^2 •^0'i.i2)(n,i2)( f )^- D (ii,i2)(fci,fc2)(*) ( 19 ) 



11,12=0 



from which the superoperator C(t) at time t can be de- 
termined. This is the main result of the paper. 

In the following we will apply this general algorithm 
into three physically interesting examples. 

Example A. - Let us consider a two- level system (a two- 
level atom, a spin- 1/2, or a qubit) with a two-dimensional 
Hilbert space Tt s spanned by two vectors |1) and |0). In 
order to specify the Liouvillian superoperator £(t) for 
the two-level atom we have to know the time evolution 
of four initial states specified by Eq. (jl5|) . Let us assume 
that from the measured data it is found that these states 
evolve as 



-rt 



e -rt ie -rt/2 



2 I -ie- Tt ' 2 2-e- rt ' ' 



1 / e~ Tt e' rt / 2 

2 { e~ rt / 2 2 - e- rt 



(20) 



Now we can apply our reconstruction scheme and we find 
for the matrix G 



G 



(jl,32)(ki,k 2 



)(*) = 



(t) the expression [[| 




(21) 



This matrix corresponds to the Liouvillian which defines 
the master equation 

d L 

-j t P = £f> = 2 [2<5-_/5(T + - v+d-p - pa+a-\ , (22) 

describing the decay of a two-level atom into a zero- 
temperature reservoir ||. The Liouvillian in Eq.(p2|) is 



time independent which reflects the fact that the state of 
the reservoir does not change in time under the influence 
of the system. 

Example B. - Here we will reconstruct the Liouvil- 
lian superoperator for the master equation describing 
the time evolution of a single two-level atom interact- 
ing with a single-mode electro-magnetic field in an ideal 
cavity. The corresponding Hamiltonian in the dipole and 
the rotating-wave approximations reads B 



H = uo, 



+ uja'a + A(fi+a + <T_a T ), 



(23) 



where A is the atom- field coupling constant. We assume 
that the atomic transition frequency (co A ) is on the res- 
onance with the filed frequency (u>). The operators a' 
and a are the usual photon creation and annihilation 
operators, respectively, with [a, a'] = 1. If the atom 
and the field are initially prepared in states \^(to)) A — 
c o |0) +ci|l), and \^{t )) F = Y,k=o e k\ k ) = \ a ), respec- 
tively, then at time t the atom-field state vector \^(t)) A _ F 
reads 

\^(t)) A _ p =co^2(cosT k \k)\0)-isinT k \k-l)\l)) 

k 

+ci5^(cos7* + i|A:}|l} -isinr fc+1 |fc+l)|0)), (24) 



where Tfc = Xy/kt. Utilizing Eq.(|l9|) we can determine the 
Liouvillian superoperator which governs the dynamics of 
the atom. Here £(t) explicitly depends on the initial 
state of the cavity field. Let us assume a particular case 
when the field has been prepared in the Fock state \M). 
With this initial state the matrix (Oh takes the form 



D 



)(*) 




(25) 



where £ = cos 2 (A£VM) and & = cos 2 (Xty/WTT) ■ The 
determinant of this matrix det[D] — £o£i(£o + £i — 1) is 
equal to zero only at discrete moments so D is invertible 
and we can utilize Eq.(|l9|) from which we find 



-7i 








7i \ 





-72/2 













-72/2 





73 








-73 / 



with the time-dependent parameters 7i(t) given as 



(26) 



71 (*) 



2 A ( \/Msin(2AtV / M) sin 2 (AtVM + 1) + y/M + 1 sin(2AiVM + 1) cos 2 (Xty/M) 



[cos(2AtVM) + cos(2A*VM + 1)] 



3 



72 W = 

73 (t) = 



(\/Msin(2Ai%/M) cos 2 (AVM + 1) + \/M + 1 sm(2Xty/M + 1) cos 2 (At\/M) 



(2XtVM) + cos(2A*VM + 1)] 



cos 2 (AtVM) cos 2 (A VAX + 1) 
2A (\/A7sin(2AtVM) cos 2 (AtVM+ 1) + VM + 1 sin(2At\/M + 1) sin 2 (AtVM) 



(27) 



From the solution (|26|) it follows that the Liouvillian 
superoperator is explicitly time-dependent which reflects 
dynamical response of the environment (i.e. the cavity 
field). The master equation (||) with C(t) specified by 
Eq.(Bq) can be written as 

7i(t) 
2 

r0) 
2 

73 (*) 



[<X + <7_/5<7_<7 + + <7_<7_|_p<7_|_<7_] 



25-4 



p(T_ — CT_(T +( — p(J — 0~ _|_J 



(28) 



with the coefficients 7j(i) Eg. (J27|) and ?y(i) = 72(0 — 
71(f ) — 73 it). One can check that pA(t) obtained from 
Eq. ( p4|) is the solution of the master equation (|2g| ) . We 
note that if the cavity field is initially in the vacuum 
state (M = 0) then the master equation (S) takes 
the form (|2^) but with the time-dependent "decay" rate 
r -> 71 (*) = 2AtanAi. 

Example C. - Finally, we consider a single two-level 
atom coupled to K modes of the electro-magnetic field 
in a one-dimensional cavity of the length L. The spec- 
trum of modes is discrete with frequencies u) k — knc/L. 
The corresponding total Hamiltonian in the dipole and 
rotating-wave approximations reads [J7J 



H 



K 

fe=i 



K 



(29) 



fc=i 



The field is assumed to be initially in the vacuum state. 
By applying our algorithm we find the master equation 
for the atom to be of the form ( j22] ) except the "de- 
cay" rate V — ► ^{t) is now explicitly time dependent, ft 
can be expressed in terms of the "measured" probability 
Pit) — (l\p A it)\l) that the upper atomic level is excited: 



lit) 



p{ty 



(30) 



In Fig.|l| we present the time evolution of P(t) and 7(t) 
obtained with the help of numerical diagonalization of 
the Hamiltonian (|2^). From our results it follows that 
7(t = 0) = but as soon as the atom starts to radiate 
the function 7(f) starts to grow and after a short time it 
takes the constant value Y = 27rA 2 d e //(w) given by the 
Golden Fermi rule |Q. At this stage the atom radiates 
exponentially and two wave packets propagating to the 
left and the right cavity mirrors are irradiated. These 
packets are reflected by mirrors at t = L/2c and they 
"kick" back the atom at t = L/c. At this point the atom 



is essentially in its ground state and the reflected waves 
packets (environment) force it to absorb energy, i.e. the 
atom does not decay exponentially anymore. This is the 
reason why during the recurrence of the atomic inversion 
7(t) rapidly changes and takes negative values. 




FIG. 1. The time evolution of the decay rate 7(4) (thin line) 
and the population of the excited atomic level P(t) (thick 
line). We assume the atom to be in the center of the 1-D 
cavity, so it is coupled only to the odd modes (i.e. \2k = 0). 
We assume L — 2ir and c = 1 so that (jJ2k+i = k + 1/2, and 
A2fc+i = A = 0.3. The effective density of modes which in- 
teract with the atom is d e //(^) = L/2cn = 1. Therefore the 
decay rate V = 27rA 2 d e// (w) ~ 0.564. We consider K = 400 
modes of the field initially in the vacuum state and the atom 
(with lu a — 101) in its upper state |1). 
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